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ABSTRACT 

A formula is derived for the rate of thermal atmospheric escape, valid, and asymptotically 
exact, at low Knudsen mimber. 



1. The thermal escape problem 

Consider a planetary atmosphere without stir- 
ring and radiative heating or cooling above the 
nominal surface. The atmosphere cannot be static 
since a static atmosphere reaches thermodynamic 
equilibrium. The Boltzmann distribution would 
then give a finite density at infinite distances, 
which is unphysical. 

The atmosphere must be in a state of perma- 
nent escape. We want to calculate the rate of es- 
cape and the resulting temperature and density 
profiles. 

The thermal escape problem has a long history 
(Jeans 1904, Parker 1958, Hunten 1982). But 
in recent applications to planetary bodies, there 
has been disagreement on which model - hydrody- 
namic or free molecular flow - should be applied 
(Johnson 2010). We believe that all issues have 
been finally settled by a direct molecular dynam- 
ics simulations of Volkov et. al. (2010). Here we 
just show how the Volkov et. al. (2010) results 
can be obtained analytically. 

It would seem that the escape rate calculation 
is only possible via molecular dynamics or Boltz- 
mann numerical simulation, because the exobase 
region, where the transition from hydrodynamic to 
free molecular fiow occurs, cannot be treated an- 
alytically. But we show, for low Knudsen number 
at the surface, that the exobase boundary condi- 
tions are either irrelevant or can be deduced due 
to large overlap between hydrodynamic and free 
molecular fiow regions. This allows one to calcu- 
late the atmospheric escape using hydrodynamics 
and to derive a formula for the escape rate. 

For clarity, we approximate the molecular in- 
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Fig. 1. — Escape rate for Kn ~ IQ^^, Cp — 5/2, 
a — \. Je < 2.5 is the Parker regime, 2.5 < Je < 
24 is the Fourier regime, Je > 24 is the Jeans 
regime. 
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teraction potential by a power law (this includes 
hard-sphere gas as a limiting case). Then, by scal- 
ing arguments, the escape rate must be given by 



Ma = Mescapc(Je,Kn). 



(1) 



Here we measure the escape rate by the radial ve- 
locity at the surface v, and the velocity is repre- 
sented by the isothermal Mach number 



Ma: 



The Jeans number is 



Je; 
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where G is the gravitational constant, M is the 
planet mass, m is the molecule mass, T is the sur- 
face temperature, R is the surface radius. The 
Knudsen number is Kn 
free path. We define 



where A is the mean 



Kn 



(4) 



where k is the thermal conductivity, p is the den- 
sity 

The escape rate formula is given in §2. The 
formula includes three branches. In §3-5, we give 
the recipes for calculating the branches. In §6, we 
explain the origin of the recipes, describe numeri- 
cal hydrodynamics which confirms the escape rate 
formula, and compare our results to the molecular 
dynamics simulations of Volkov et. al. (2010). 

2. The escape rate formula 

In the low-Knudsen limit. 



Mparkor(Je), 
Mcscapo = { MFouricr(Je,Kn), 
Mjoans(Je), 



Cp < < ln(Kn~ 
Je > ln(Kn-^). 

(5) 

The formula is 



Here Cp is the heat capacity, 
asymptotically exact (Appendix). 

Mparkcr is just the Parker wind (Parker 1958) 
escape rate. It depends on Je and Cp. For example, 
for Cp = 5/2, 



Formula ([5]) limits the Je range over which the 
Parker rate applies. The Parker regime is de- 
scribed in §3. 

Mjcans is the Jeans (1904) escape rate. It de- 
pends on Je only: 

(Je + 1) exp(— Je) 
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Formula ([5]) extends the Kn range over which the 
Jeans rate applies. The Jeans regime is described 
in §4. 

Mpouricr is proportional to Kn, with the Je- 
dependent coefficient of proportionality: 



MFouricr(Je,Kn) =F(Je)Kn. 
Here F(Je) depends on Je, Cp, and 



d\nT 



For example, for a = 1, 
F(Je) = 



Je 



(8) 



(9) 



(10) 
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For uniformity of notation, we call it the Fourier 
regime. The Fourier regime is described in §5. 

To illustrate the formula. Fig. 1 shows the es- 
cape rate as a function of Je, for Kn = 10~*, 
a = 1, and Cp — 5/2. The different regimes were 
extended to the intersection points. Fig.l is in 
qualitative agreement with molecular dynamics re- 
sults of Volkov et. al. (2010). Quantitative com- 
parison is done in §6, where we also show that 
a simplified dissipative hydrodynamics model re- 
produces the molecular dynamics results to a good 
accuracy. 

In the next three sections we calculate the es- 
cape rate postulating the properties of the flow. 
The assumptions made can be shown to be self- 
consistent. But, in fact, what makes the results 
trustworthy is that they do match the numerical 
hydrodynamics results described in §6, and they 
also match the molecular dynamics simulations of 
Volkov et al (2010). 

The Parker regime is a high- velocity ideal flow 
with negligible thermal conductivity. The Fourier 
regime is a flow with algebraically small velocity 
and with thermal conduction balancing the ideal- 
hydrodynamic energy flux. The Jeans regime is a 
flow with exponentially small velocity, such that 
thermal conduction dominates. 
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3. Parker regime 

For Je < Cp, Kn <C 1 one can use ideal hydro- 
dynamics 0. To write equations in the simplest 
form, we use scalings to set M = m = R = 1, and 
also p — 1, T ~ 1, a.t the surface R = 1. Now 
we need to calculate the outflow velocity at the 
surface vq — Mparkor, for any Je = G < Cp. 

We write stationary ideal hydrodynamics as 
flux conservation, energy flux conservation, and 
adiabaticity 

pvr'^ = vo, 

= 4 + Cp~G, (11) 

From eas. pT|) 

— + cp(— = — +Cp-G+— . 12 

This "algebraic" equation for v must have solu- 
tions for all r > 1, such that the corresponding 
density p = vn/{r^v) goes to zero at large r. This 
requirement limits the range of possible values of 

Vq. 

For Cp < 3, and also for Cp > 3, G < 
2cp/{cp — 1), the smallest possible value is vq = 
\/cp/{cp — 1), in which case the sonic point is 
at the surface. For Cp > 3 and larger G, 
G > 2cp/{cp — 1), the sonic point moves out to 
r > 1, and vq decreases. 

It is the smallest possible vq that we call 
Mparkcr- Timc-dcpendent ideal hydrodynamics 
simulations and simulations of Volkov et al (2010) 
confirm that the flow does approach the smallest 
Vo allowed by stationary ideal hydrodynamics. 

The above results can be written as 

Mparkcr- ./^^, (13) 

for Cp < 3, and also for Cp > 3, Je < 2cp/ (cp — 1). 
For Cp > 3, 2cp/{cp - 1) < Je < Cp, Mparkcr de- 
creases, reaching zero at Je — Cp. The precise 
value of Mparkcr IS then given implicitly by an 



More precisely, dissipative effects must become negligible 
after a possible boundary layer is passed. Whether the 
boundary layer is present, depends on how the atmosphere 
is fed. If the boundary layer is present, we simply move the 
nominal surface up, beyond the boundary layer. 



"algebraic" equation which can be derived from 
eq.dlllEI. 

4. Jeans regime 

Fix the Knudsen number and let the Jeans 
number increase. Then the escape rate approaches 

M,.an. - (14) 

V 27r 

- the standard Jeans rate - density normal- 
ized flux of positive-energy positive-radial- velocity 
molecules, calculated for Maxwellian distribution. 

Ea. (ll4p is an exponentially accurate approxi- 
mation (9% error at Je = 2.5, 0.008% error at 
Je — 10) of the true escape rate for collision- 
less molecular outflow with Maxwellian injection 
at the surface, which in our notations is 

f„°° udu f°° dvfv 

TV J- JO J — OO ■> rN 

iVl Jeans — /.oo , noo , „ • 

Jo 'uduj_^dvf 

Here v is the radial velocity, u is the absolute value 
of the tangential velocity, / cx 9 exp{~{u'^ + v'^)/2) 
is the truncated Maxwellian. 9 = ii simultane- 
ously the energy is positive, (m^ -|- w^)/2 — Je > 0, 
and the radial velocity v is negative; otherwise 
9 = 1. 

5. Fourier regime 

Fix the Jeans number Je > Cp and let the Knud- 
sen number decrease. The escape rate then ap- 
proaches 

Mpouricr = F(Je)Kn. (16) 
The coefficient i^(Je) is calculated below. 

In the Fourier regime, the negative hydrody- 
namic energy flux carried by the flow is almost 
exactly compensated by the thermal conductivity: 

r^KT' = ^{c„T--). (17) 
r 

Here $ = pvr^ =const is the flux per steradian, 
M = m = 1 was set by scaling, prime denotes 
the r-derivative, kinetic energy and viscosity are 
negligibly small. 

Using scalings, we put T = 1 at r = 1, and then 
solve eq. pT]) for r > 1. We adjust $ in such a way 
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as to get the temperature going to zero at large 
radii (otherwise the flow either stahs or becomes 
singular at finite radii, never reaching an exobase, 
§6). We want to solve the problem for arbitrary 
a = dlnK/dlnTE- 

First consider just the case a — I. One can 
check numerically or analytically, that there is 
only one solution of eq. (jl7|) with correct asymp- 
totic behavior at large r: 



r 

Substitution into ca. (|T7| then gives 
^ = (G-c,)<i>. 

Or, in difi^erent notations, 
F{Je) 



(18) 



(19) 



(20) 



For arbitrary a < 1, numerical integration of 
ea. ((T7)) shows that there is still a unique choice of 
$ which gives a non-singular zero-at-infinity tem- 
perature. One can show that for < Je — Cp <C 1, 



F 



1 



Je 



while for Je » 1, 



1 -I- a Je 

Joining the asymptotics as follows 
F(Je) ^ 



/(Je) 



1 



(Je - Cp)/(Je) 

1 — a 3e — Cp 
2 Je - 0.55c„ 



(21) 



(22) 



(23) 



(24) 



gives an error less than 0.3% for all relevant a, Cp 
and Je. 

6. Hydrodynamics 

To show that the escape rate formula ([5]) is 
asymptotically exact at Kn — > 0, and to improve 
the accuracy in the non- asymptotic regime, we 



•^a = 1 is a good approximation for real gases at relevant 
temperatures; a = 0.5 for the hard-sphere gas, and we 
want to reproduce the results of Volkov et. al. (2010) 



use stationary non-ideal hydrodynamics. The flux 
conservation, the energy flux conservation, and the 
momentum equations are 

Q = H^+ CpT - f ) ™ l^7^r^v{v' - f ) - r^nV, 

pvv' = -{pry + l^ir^viv' - f ))'• 

(25) 

Here Q are constant fluxes, 77 is the viscositjQ. 
The boundary conditions at r = 1 are 

/9=l,T = l,w = $ (26) 

Given G = Je and k(T = 1) = Kn, we must find 
a pair ($,(5) that gives the "true" flowE]. Which 
flows are "true" is decided by the boundary con- 
ditions at the exobase, and it would appear that 
one cannot proceed without a molecular dynamics 
simulation of the exobase region. 

But it turns out that for Kn ^ 1 the mere 
requirement that the exobase exists allows to cal- 
culate the escape flow in the Parker and Fourier 
regimes. In the Jeans regime, this trick no longer 
works - there are many outflows with an exobase. 
What allows to solve the problem in the Jeans 
regime is a large overlap between hydrodynamic 
and free molecular flow regions. 

Parker regime - Take G < Cp, and, as always, 
K ~ Kn ^1. If the initial velocity is large, 
vq = ^ ^ 1, the non-ideal hydrodynamics (l^5t de- 
generates into the ideal hydrodynamics PT|) . and 
we get the Parker escape. 

Only for small velocities, vo ^ k, dissipative 
effects can be significant. But numerical integra- 
tion of ([^5]) shows that for small vq the exobase 
is never achieved. Kn stays small at all radii, no 
matter what energy flux Q one chooses. 

What sets the boundary of the Parker regime 
to Je = Cp is the sign of the ideal hydrodynamic 
energy flux in the small-velocity limit 



Qh = *(cpr--). 

r 



(27) 



* We have assumed that the second viscosity vanishes. This 
is true for the hard-sphere gas and allows precise compar- 
ison to the results of Volkov et. al. For polyatomic gases, 
the second viscosity must be included, but it turns out that 
the viscosity effects are negligible at low Knudsen, unless 
Je is very close to Cp. 

^ To solve I I25II we also need v'{l). But this boundary con- 
dition is irrelevant. Wrong ^^'(l) would result in a bound- 
ary layer, and we have agreed to place the nominal surface 
above the boundary layer. 
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For Je < Cp this is positive. Now we show that 
this means no exobase for small vq. 

For small vq , the density is approximately given 
by the equilibrium condition 

- (PT)' - ^ = 0. (28) 

At moderate G, the only way to reduce the density 
enough to get an exobase, is by decreasing the 
temperature at large distances. If the temperature 
remains finite at infinity, the density saturates as 

p^exp(^). (29) 

At small Vq , we can drop the terms from the 
energy flux expression in (1^51) . and get 

r^nT' ^<i>(cr,T--)-Q. (30) 
r 

For the temperature to decrease, we must take a 
positive fiux Q. But finite Q then dominates the 
r.h.s. of ca. (l30l) at large radii, and we get 

r^nT' = -Q. (31) 

This gives 

Tcxr^TT^. (32) 

This temperature decreases - but too slowly. For 
a > 0, eas. ipSl give a density growing at in- 
finity. 

In summary, small-velocity fiows have no 
exobase, and large-velocity fiows are ideal, giv- 
ing the Parker fiow. The lack of collisions beyond 
the exobase cannot change the escape rate of the 
highly supersonic outfiow. 

Fourier regime - At fixed Je > Cp, as Kn de- 
creases, the set of pairs ($, Q) which give solu- 
tions with an exobase shrinks to a single point 
(MFouricr,0). If ($, Q) ^ (Mpouricr , 0) , the fiow ei- 
ther becomes singular at a finite radius or goes to 
infinite radii, without ever reaching local Kn ~ 1. 
So in this case too, we can claim that Mpouricr is 
an exact asymptotic even without knowing what 
happens at the exobase. 

The only way to reach an exobase is to set Q — 
in ea. ([30| . and then tune <I>, as described in §5. 

Jeans regime - Fix Kn and let Je increase. 
The outflow velocity becomes so small that ther- 
mal conduction in eq. (j30p dominates and we get a 



nearly constant temperature profile. Eq. ([28)) then 
gives an exponential density. These temperature 
and density are also predicted by the free molec- 
ular outflow at high Je. We see that the two de- 
scriptions agree. 

We can therefore use hydrodynamic descrip- 
tion well above the point r = where the 
density scale height becomes comparable to the 
mean free path, Kn(i?e) — Je~^(i?e)- But for 
Kn(r) ^ Je~^(r), we can also use the free molec- 
ular outfiow results to calculate the fluxes ($, Q). 
Now, we use hydrodynamics to integrate back to 
r = 1. One can show, or just check numerically, 
that this procedure, for Je ^ InKn^^, gives the 
Jeans rate (to leading order in Je in the prefactor). 

Simplified hydrodynamics - Both the Fourier 
and the Jeans regimes can be described simulta- 
neously, in a much simplified version of hydrody- 
namics which we give below. One can then hope 
that even between the asymptotics this hydrody- 
namic description will work (it does work for hard 
spheres, see below). 

The model has one adjustable parameter Kuq ^ 
1. We set T = 1, p = 1, Kn = Kng at r = 1, and 
use the free molecular flow to calculate $ and Q 
at r = 1: 

^_ e-g(G + l) e-g((cp-|)G + Cp-l) 

\/2^ ' V2^ 

(33) 

We then integrate eg. ((30)) backwards in radius, to 
r < 1. Simultaneously we integrate eq.(|28l) and 
calculate Kn(r). The integration is terminated 
when the desired value of Kn is reached. We ad- 
just G so that the desired value of Je is reached 
simultaneously with the desired value of Kn. This 
gives the escape rate for these Kn and Je. One 
also gets correct temperature and density profiles 
for G ^ 1 and also for r ^ 1 at any G. We have 
checked, by numerical integration of the full sys- 
tem ([25)) . that neglecting viscosity and inertia is 
justified, again for G 2> 1 and also for r -C 1 at 
any G. 

We have used the simplified hydrodynamics 
(I30I33P to establish the logarithmic accuracy of 
the escape rate formula ([5]) in the Fourier- Jeans 
transition region (Appendix). 

Comparison with Volkov et. al. (2010) - Now 
compare the predictions of our escape rate formula 
eq.(IS|) to the molecular dynamics simulations of 
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Volkov et. al. (2010). We put Cp = 5/2, a = 1/2. 
We also remember that our definition of Kn is dif- 
ferent by a constant factor from the conventional 
hard-sphere Knudsen Knhs — 0.34Kn. 

We find that for Je = 3,10,15, Knhs = 
10~^, 10~^, 10"'' the simplified hydrodynamics 
P0|33|) (with Kno = 3) predicts the escape rates 
which agree with the resuhs of Volkov et. al. 
(2010) to about 10% for all points but one. 

Let us discuss this discrepancy. For Je = 3, 
the escape rate formula dU gives for Kuhs = 
10-2,10-3,10"'' (using eqs.dlllESllHl)) 



Ma 
Kuhs 



= 6.4, 6.4, 6.4 



The simplified hydrodynamics p3|) gives 



Ma 
Knhs 



= 5.6, 6.0, 6.2 



Fig.4 of Volkov et al (2010) gives 



Ma 
Knhs 



= 5.0, 6.8, 4.1 



(34) 



(35) 



(36) 



The simplified hydrodynamics works as expected 
- as we go deeper into the Fourier regime, pre- 
dicted Ma approaches Mpourior. The first two val- 
ues of Volkov et. al. are also reasonable - the 
error of Mpouricr decreases as we get deeper into 
the Fourier regime. But the last point should have 
been even closer to the Fourier value of 6.4. Some 
50% error is unexpected. 

We believe that it is the molecular dynamics 
results at Je — 3, Knhs — 10^^ which have a large 
error. We have checked, by integrating the full sys- 
tem (|25p numerically, that a hydrodynamic fiow, 
starting at Kuhs = 10""*, Ma = 4.1 x lO"", Je = 3, 
never goes above Kuhs = 1.5 x 10^"^ for any Q. It 
seems unlikely that molecular dynamics can no- 
ticeably deviate from hydrodynamics at such low 
Knudsen numbers. 

7. Summary 

Formula ([5]) gives the thermal escape rate at 
low Knudsen. 

The escape occurs either in Parker, in Fourier, 
or in the Jeans regime. In the Parker regime, 
we get a supersonic ideal hydrodynamical wind. 
The temperature and density profiles are given by 



the Parker wind model (§3). Both temperature 
and density decrease algebraically. In the Jeans 
regime, we have a nearly isothermal slowly escap- 
ing atmosphere, with exponential density. In the 
Fourier regime, thermal conductivity balances the 
ideal energy flow (§5). The temperature and the 
density decrease algebraically. 

I thank Robert Johnson for useful discussions. 

Appendix 

The escape rate formula (5) is asymptotically 
exact in the following sense: 



Je = const < Cp, Kn 



M, 



Mparkcr 



1, 



(37) 



Je = const > Cp, Kn — ^ 



Mc 



-Mpouricr 



1, (38) 



Kn = const, Je — >■ cxd 



Je > Cp, Kn 



M, 



escape 



Mj 

cans 
In(MoscapG) 



1, 



ln(min(MFourior, Mjcans)) 



(39) 



1. 



(40) 

In words: (i) As Knudsen decreases at fixed 
Jeans, the Parker and Fourier formulas approach 
the true escape rate everywhere but near the 
Parker-Fourier transition. As Knudsen decreases, 
the Jeans number width of the Parker-Fourier 
transition region shrinks to zero, (ii) As Jeans 
increases at fixed Knudsen, the Jeans formula ap- 
proaches the true escape rate, (iii) In the Fourier- 
Jeans transition interval, both formulas are valid, 
but only with logarithmic accuracy. 
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